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In the mixed state of a d-wave superconductor, Bogoliubov quasiparticles are scattered from magnetic vortices via a 
combination of two effects: Aharonov-Bohm scattering due to the Berry phase acquired by a quasiparticle upon circling 
a vortex, and effective potential scattering due to the superflow swirling about the vortices. In this paper, we consider 
the Berry phase contribution in the absence of superflow, which results in branch cuts between neighboring vortices 
across which the quasiparticle wave function changes sign. Here, the simplest problem that captures the physics 
is that of scattering from a single finite branch cut that stretches between two vortices. Elliptical coordinates are 
natural for this two-center problem, and we proceed by separating the massless Dirac equation in elliptical coordinates. 
The separated equations take the form of the Whittaker-Hill equations, which we solve to obtain radial and angular 
eigenfunctions. With these eigenfunctions in hand, we construct the scattering cross section via partial wave analysis. 
We discuss the scattering effect of Berry phase in the absence of superflow, having considered the superflow effect in 
the absence of Berry phase in a separate paper. We also provide qualitative comparison of transport cross sections 
for the Berry phase and the Superflow effects. The important issue of interference between the two effects is left to 
future work. 



I. INTRODUCTION 



Massless Dirac quasiparticles are the low energy excitations of a d-wave superconductor. In the vortex state, low 
temperature transport properties, such as longitudinal thermal conductivity and thermal Hall conductivity, can be 
explained by studying the scattering of these quasiparticles from magnetic vortices. Quasiparticles scatter from 
vortices via two basic mechanisms: a circulating superflow and a Berry phase factor of (-1) acquired upon circling a 
vortex. Scattering due to superflow current in a single vortex (without the Berry phase effect) has been calculated in 
Refs. Q] and El Although the transport properties calculated in Ref. 1 capture the essential physics qualitatively, the 
Berry phase effect becomes important for the higher field (lower temperature) regime where the deBroglie wavelength 
is comparable to the distance between vortices. It is therefore of interest to consider the effect of Berry phase on the 
quasiparticle scattering. The influence of the Aharonov-Bohm effect due to an isolated vortex line on the quasiparticle 
states has been studied in the Refs [3] and |U Ref [3] also obtains the quasiparticle wavefunction and density of states 
at distances from the vortex large compared to the pentration depth (r >> A). In this paper, we consider (r << A) 
and work in the limit A — > oo, such that the magnetic field is constant across the sample. This limit permits 
the application of a singular gauge transformation that encodes the Berry phase effect in the form of anti-periodic 
boundary conditions on the quasiparticle wave function. In contrast to conventional superconductors, the density of 
states at low energies in d-wave systems is dominated by contributions which come from the regions far from the coreiP 
and are associated with extended quasiparticle states with momenta close to the nodal directions. This allows us to 
neglect the vortex core physics by setting the coherence length (£ — > 0). This conclusion, based on the semiclassical 
approacrPEl, has been confirmed by numerical analysis-^ of the BdG equations for a single isolated vortex line (in the 
limit A — > oo). Note that the calculations presented in|7|also point to the absence of truly localized core states or 
any resonant levels in the pure d-wave case, though such states were observed in numerical simulations^. AB effect 
on quasiparticle excitations in macroscopic superconducting rings has been studied in detail in Refs [9l and 1101 In this 
work, we calculate the scattering contribution due to the Berry phase effect. We take the following path. As discussed 
in Ref. [5] and summarized in Appendix [X] we apply a singular gauge transformation to the Bogoliubov-de Gennes 
equation and shift the origin of momentum space to the location of one of the gap nodes. This reduces the problem 
to that of an (anisotropic) Dirac fermion scattering from an effective non-central potential (due to the superflow) in 
the presence of antiperiodic boundary conditions (a consequence of our gauge choice) . For simplicity, we neglect the 
anisotropy of the Dirac dispersion by taking Vf = v%. Since we are only interested in the Berry phase contribution, 
we neglect the presence of the effective non-central potential, which further reduces the problem to that of a massless 
Dirac fermion scattering due to the antiperiodic boundary conditions. Within the single vortex approximation, the 
antiperiodic boundary condition manifests as a semi-infinite branch cut terminating at the vortex core. In Sec. |H| we 
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study the single vortex scattering of quasiparticles due to this semi-infinite branch cut (without superflow current) 
and obtain a divergent differential cross section in the forward direction. The divergent nature of this cross section is 
unphysical and requires that we regularize this semi-infinite branch cut. In real situations these branch cuts terminate 
on the cores of the neighboring vortices and are finite in nature. Thus, considering a pair of vortices as our scatterer 
captures the finite branch cut and regularizes the Berry phase effect. Elliptical coordinates serve as a natural choice 
for this problem with two vortex cores. The presence of a finite branch cut between the two vortices manifests itself 
as a boundary condition on the wave function spinor across the line segment joining the two foci of the ellipse. In 
elliptical coordinates (/i, v), we simply write 

U=o=(-i) B U=o (1) 

We can impose the Berry phase condition using parameter B (B = 0,1). When B=0, there is no branch cut 
between the vortex cores. We can turn on the Berry phase (branch cut) between the vortex cores by setting B=l. 
Thus, neglecting the anisotropy of the gap nodes, our problem reduces to that of the scattering of massless Dirac 
quasiparticles due to this finite branch cut. A similar setup for the scattering of quasiparticles due to the Berry phase 
has been considered by Melikyan and Tesanovic in Ref. 1111 Their approach was to construct scattering solutions to 
the 2d Dirac equation in elliptical coordinates from solutions to the Klein-Gordon equation (separated in elliptical 
coordinates) using self-adjoint extensions. But in going from Klein-Gordon to Dirac, the solutions are not separated 
in elliptical coordinates anymore. This results in not being able to construct all the self-adjoint extensions. In this 
paper, we avoid the issue of self-adjoi nt ext ensions by taking a different approach. We dire ctly separate the (2+l)d 
Dirac equation in elliptical coordinate d 12 * 13 !. The separation of variables is done in Sec. IV The separated equations 



are Whittaker Hill equations (WHE^i 2 ^! In Sec. [v]we solve the separated equations to obtain eigenfunctions^^U. 
In Sec. VI we develop an expansion for the incident plane wave (representing quasiparticle current) in terms of the 
separated solutions of WHE. In Sec. | VII| we construct the scattering amplitude (as a sum of partial waves) from the 
asymptotic form of the scattered wave. We impose the boundary condition Eq. on the full wavefunction spinor 
and calculate the phase shifts for each partial wave. We show that for B=0 (no branch cut) there is no quasiparticle 



scattering. In Sec IX we turn on the Berry phase by setting B=l (with branch cut between cores) and obtain a 
non-zero scattering cross section. Results and analysis of the quasiparticle scattering cross section in the presence of 
Berry phase is presented in Sec. [Xl Conclusions are discussed in Sec. |XI| 



II. BERRY PHASE SCATTERING OF INCIDENT PLANE WAVE IN SINGLE VORTEX 
APPROXIMATION (WITHOUT SUPERFLOW) 



Quasiparticles scatter from vortices via both the circulating superflow and the Berry phase factor of (-1) acquired 
upon circling the vortex. This phase is encoded in the antiperiodic boundary conditions imposed on quasiparticles in 
our chosen gauge (see Appendix [A] for details). In this section, we neglect the superflow by setting P s = and consider 
only the Berry phase contribution. That is, we consider the scattering of quasiparticles due only to the presence of 
antiperiodic boundary conditions. Furthermore, we neglect the anisotropy of the Dirac dispersion and take Vf = V2 
(a = 1). As we shall see, the antiperiodic boundary conditions yield the Aharonov-Bohm interference effect of an 
enclosed 7r-flux. 

We consider the isotropic Dirac Hamiltonian 



H = Vf [r 3 p x + npy] 



(2) 



and seek solutions to the Bogoliubov-de Gennes equation, = E<&. We express the quasiparticle wave function 
as a linear combination of angular momentum eigenstates which satisfy = Since we have neglected the 

superflow, there is no effective potential and the general solution is easily found to be 



1 



■i + B ^ 



1 



(3) 



where and are complex constants and </ M ±i and are Bessel functions of argument p = kr. However, 

rather than imposing periodic boundary conditions by requiring that p = n + 1/2 with n = integer, here we shall 
impose antiperiodic boundary conditions, 



$(r, (j) + 2n) = -$(r,(/>) 



(4) 
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$,(r) = e^e** ( ^ | ) 7 = ±1 (5) 



by requiring that /i = integer. The radial functions are therefore half-integer Bessel functions rather than integer 
Bessel functions. The coefficients, and B^, are determined by satisfying boundary conditions at both long and 
short distances. 

At long distances, we require an asymptotic wave function that is equal to the sum of an incident plane wave, <!>i, and 
an outgoing radial wave, $ s . In terms of the current functional discussed in Appendix |a| j[$] = Vf^^x + riy)$, 
we require that ~ k and j[$ s ] ~ r. In the presence of antiperiodic boundary conditions, we seek an incident 
wave of the form 

and a scattered wave of the form 

where k = (k, 9), r = (r, cf>), and ip = <fi — 9. In analogy with the problem of Aharonov-Bohm 2 ^ scattering from an 
enclosed magnetic flux, we can say that 7 = — 1 corresponds to an effective 7r-flux while 7 = +1 corresponds to an 
effective (— 7r)-fmx. Since these two cases are equivalent, the choice of 7 = ±1 is arbitrary. In the asymptotic limit, 
the half-integer Bessel functions take the form 

where rj^ — 1 for fi > and rj^ — (— 1)** for \x < 0. Proceeding along the lines of Ref. |5J we can plug these asymptotic 
expressions into Eq. (|3| , reorganize terms, and thereby obtain a suggestive (yet still general) form for the quasiparticle 
wave function. Doing so, we find that 

^ = e nf,ik.r /cos I A . „ H e lkr ( cos t \ . , ,e~ lkr 

where 



s Vs E »(*>) = E w 

and a M and are complex constants defined via 

Afj, - i-B^ = i^e"*"* (e-^/2 + 6„) (10) 



A M + iB^ = -T^e-*"' (e^/2 + (-l)"o M ) . (11) 

If the plane wave is to be the only incident wave, we must eliminate the incident radial wave by requiring that a M = 
for all /.i. With this restriction, and 6 M are related to via 

B M = i(A,+* 7 A°) 6|t = e -»* - (12) 
where A° = i^-^e'^ 6 /2. The asymptotic wave function then takes the desired form 

* = e^*(™j) + f(^( CO i) (13) 
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and the differential cross section is given by 



(14) 



We can now determine the by imposing appropriate boundary conditions at the origin. As discussed in Ref. [3| the 
most restrictive condition is that the current through the origin (a point of zero area) must be zero. More precisely, 
we consider a semicircle of radius e, oriented about the 9 direction, and require that the total current passing through 
it, Ig, vanishes as e — >• 0. If $(/? — >• 0) ~ p a , then Ig ~ e 2a+1 . Thus, to ensure that Ig does not diverge at the origin, we 
must eliminate all terms in Eq. (|3b which diverge faster than p -1 / 2 as p — >• 0. Since the half-integer Bessel functions 

" M=F 5, this clearly requires that 



exhibit the small-argument behavior, J^±k ~ p M± 2 and Y^±i 
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= for p > Af, = for fi < 

The condition for p = is more subtle. Enforcing the above, the resulting wave function is dominated, as p 
the terms which diverge exactly as p^ 1 / 2 . We therefore find that 



*(/> -> 0) 



2 

7rp 



(A - iB ) 



| + ^( J 4 + iB ) 



(15) 
0, by 

(16) 



From Eq. (12 1, we know that A + iB = —vyA^ where Aq = e !7r / 4 /2. Furthermore, we can define a complex 
constant, j3, such that Aq — iBo = — {jAq/3. With this definition, the current density near the origin takes the form 



j(p->0) = 



v f 
2np 



l)f + 2Im\0\4> 



Explicitly computing the current through the origin, we find that 



f.e+Tr/2 

Ig = lim / ed(pj(e) ■ 9 = vj 



0-7T/2 



2tt 



1 f 



■2Im[p]$' 



(17) 



(18) 



which must be set to zero for all directions 9. This requires that j3 = ±1. 
Putting everything together yields the values of our original coefficients 



A„ 




B„ 




(19) 



where 7 = ±1 and f3 = ±1. The Zi ambiguity in 7 and /3 is a consequence of the equivalence of a 7r-flux with a 
(— 7r)-flux, which cannot affect observable quantities. For /3 = 7 = ±1, 



J_ f i 7^ > 
72 I -1 < 



Plugging this into Eq. (14 1 and summing over p (with a convergence factor e l' 1 l 0+ ), yields 



da _ 1 

d<£ 27rfc sin 2 (c/?/2) 



,2%) 



(20) 



(21) 



The same result is obtained for /3 = —7. Note, however, that the above is only valid for cp ^ 0. As discussed (for the 
electron scattering case) in the original paper by Aharonov and BohnPU, as well as in an excellent review by Olariu 
and PopesciP^, our asymptotic approximations are only valid away from the forward direction. Thus, the second 
term above, which is only nonzero for (p — and is an artifact of our casual treatment of the forward direction, can 
be dropped. (For a detailed treatment of the Aharonov-Bohm scattering of an electron in the forward direction, see 
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the paper by StelitancP^ Our differential cross section therefore takes the form 



dip 2nk sin 2 ((/j/2) 

which is exactly the differential cross section for the Aharonov-Bohm scattering of an electron from an enclosed 7r-flux. 
As expected, this result is independent of our choice of 7 = ±1 and j3 = ±1 and is the same for quasiparticles about 
any of the four gap nodes. Due to the infinite range of the Berry phase effect, the total cross section diverges. However, 
the transport cross section is finite and given by = 1/nk. Since left-right symmetry is not broken in the absence of 
a superflow, the skew cross section is zero. In the zero-superflow case considered above, it was easy enough to neglect 
the subtleties associated with forward scattering in the presence of antiperiodic boundary conditions. However, if we 
were to consider the superflow and the Berry phase effects together, it would be necessary to treat such nuances more 
carefully. The first step towards that is to regularize the calculation of the cross section due to the Berry phase effect. 
This is the goal of the remainder of this paper. 



III. REGULARIZATION OF BERRY PHASE IN DOUBLE VORTEX SETUP 



Branch Cut 




Figure 1: Single vortex with semi-infinite branch cut and double vortex with finite branch cut due to the Berry phase 



The infinite range of the Berry phase effect is solely due to the isolated treatment of the single vortex (see Fig. [TJ. 
In reality these vortices are not isolated and the Berry phase effect terminates at the neighboring vortices thereby 
resulting in a finite branch cut. The simplest object containing a finite branch cut is a pair of vortices separated by 
some distance as shown in Fig. [T] Elliptical coordinates is a natural setting for studying two-center problems. We 
define elliptical coordinates via 

x = R cosh /i cos v (23) 
y = R sinh fx sin v (24) 

The presence of a finite branch cut between the two vortex cores can be imposed as a boundary condition at fi = 
(line separating foci) on the full wavefunction spinor in the following way 

-v) U=o= IK/*, f ) U=o (25) 

Parameter B is defined to take values or 1. For B=0, there is no branch cut between the pair of vortices (or foci) 
and the condition reduces to a trivial continuity condition at fi — of the elliptical coordinate system. Setting B=l, 
we turn on the Berry phase effect via non-trivial boundary condition. The presence of the branch cut captures the 
fact that when quasiparticle passes between the vortex cores it acquires a phase factor of (-1), which on its own can 
scatter quasiparticles. Our task for the rest of this paper is to capture the scattering of quasiparticles due to a finite 
branch cut between two vortices. We will follow the same prescription that we applied to the single vortex case to 
calculate the scattering cross section. The first step is to separate the (2+l)d Dirac equation in elliptical coordinates. 
The second step is to solve the separated equations to get the eigenfunctions for the separation parameter. Third 
step is to construct incoming plane wave and outgoing scattered wave in terms of phase shifts from the separated 
eigenfunctions. Fourth and final step is to impose the boundary condition for the branch cut on the full wavefunction 
spinor and calculate the scattering amplitude and differential scattering cross section. 
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IV. SEPARATION OF DIRAC EQUATION IN ELLIPTICAL COORDINATES 



In this section we separate the Dirac equation in elliptical coordinates!- 2 . Detailed calculation for separation of 
variables can be seen in Appendix [B] Dirac equation in (2+l)d is given as 



[7°5 t +7 1 ^+7 2 9,]^ = 



(26) 



We define the Dirac matrix representation to be 7 = T2, 7 1 = ir\ and 7 2 = — ir^ , where r/s are the Pauli matrices 
We define elliptical coordinates via 



x = f(p, v ) = R cosh/i cosv 
y = g(ji, v) = R sinh /i sin v 



(27) 
(28) 



Here fi is the radial coordinate and v is the angular coordinate. Here, we directly define the transformations that 
separate the Dirac equation in elliptical coordinates and leave the rigorous details to Appendix IS] 



h h 

The product of transformation matrices S and T can be explicitly written out as 

„j, _ 1 / cos I cosh I — sin I sinh ^ 

~ (cosh fx + cos v) \ sin | sinh | cos | cosh ^ 



(29) 



(30) 



ifj = STY 

Y is the transformed wavefunction spinor and is given by 



Y 



Each component of Y satisfies separable second order equation. We define a separation constant A 2 such that 

(<9 2 - ikR cos v + k 2 R 2 sir? v + \ 2 )A{v) = 
(pi +ikRcosv + k 2 R 2 sin 2 v + \ 2 )B{v) = 
(<9 2 -ifc J Rcosh^ + fc 2 i? 2 sinh 2 //- A 2 )a(^) = 



(<9 2 + ikRcosh/j, + k 2 R 2 sinh 2 fi - X 2 )/3(n) 







(31) 



(32) 



(33) 
(34) 
(35) 
(36) 



which is equivalent to 4 coupled first order equations that connect the upper and lower components of the wave 
function spinor 



{d v - ikRsinv)A(v) = i\B(u) 

{d v + ikRsvav)B{v) = i\A(v) 

{8^ — ikRsinh fj,)a(fj.) — \/3((i) 

{d^ + «fci?sinh/i)/3(/i) = Aa(/i) 



(37) 
(38) 
(39) 
(40) 



The separated equations Eq. (33), (34), (35), and (p6| are known as the radial and angular Whittaker Hill equations 



(WHE) and the upper and lower components of the spinor are connected via first order coupled equations. 



V. SOLUTIONS TO THE WHITTAKKER HILL EQUATION (WHE) 



We transform the radial and angular WHE by using the following functional transform. 

A{v) = y x e- lkR msv 



(41) 
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where y\ satisfies the differential equation known as the Ince equation^ * 17 * 22 ! . 

y 1 + 2ikR sin v y 1 + A 2 yi =0 
We do similar transformations for the radial WHE 

a{p) = x 2 e lkR cosh " 

to obtain 

x 2 + 2ifci?sinh px 2 — \ 2 x 2 = 



(42) 
(43) 
(44) 



B(v) and (3(p) are calculated from the first order coupled equations. Now we will try to solve all the Ince equations 
as an eigenvalue proble m u sing a m atrix method. We write a general form of the Ince equation from which we can 
deduce Ince equations ( 42 1 and ( 44 1 and try to find the recursions for the general equation^. 



d 2 ib dib 

-— £ + 2i uj sin Q-j- + (A 2 + uj(p + i) cos 9)ip = 
aa z do 



(45) 



if) = j/i, 6 = v, p = —i and oj = kR will yield Eq. (42 I. 



■0 = x 2 , 9 = ip, p = —i and uj = —kR will yield Eq. (44 1. 



A. Solutions to angular WHE (Matrix method) 



As often happens in the solution of differential equations with periodic coefficients, the solutions fall into four classes 
corresponding to the four types of fourier series. They may be even or odd functions of 8 and may have 2ir as their 



period or antiperiod. The four possible solutions to Eq. ( 45 1 are 

oo oo 
r=0 r=0 

= a ™ cos(r + -)9, Vm(^) = E h ™ sm(r + (47) 



r=0 r=0 



We know that the full wave function has to have 2n periodicity in v. In order to have this we have to choose 
the solutions of the Ince equations to be 2-n antiperiodic in v. This is so because the transformation matrix (ST) 
multiplying the wave function spinor is 2-k antiperiodic in v (see Appendix|B]for details of the transformation matrix), 
hence the product of antiperiodic spinor and antiperiodic transformation matrix will yield 2-k periodicity in v. Hence 



the possible solutions are Eq. ( 47 1 . Now we substitute these solutions in Eq. ( |45| and get recursion relations for the 
coefficients. Reducing WHE to Ince equations has the advantage that the Ince equation has only three term recursion 
relations as opposed to the five term recursions for WHE. The eigenvalue equation for the even solution is 



2 1 UJ P\ 10 



a Q (\e'-- + ^) + a 1 -(p-2 l ) = (r = 0) (48) 
|(p + 2riK-i + (Ae 2 -(r+i) 2 K + ^(p - 2(r + l)i)a r+1 - (r > 1) (49) 



The eigenvalue equation for the odd solution is 

6o(A 2 -i-^)+6 1 |(p-2z) = (r = 0) (50) 

U ^{p + 2ri)b r _ 1 + {\o 2 -{r+ l -f)b r + ^{p-2{r + l)i)b r+l - (r > 1) (51) 

We can see that the coefficients of ao and bo are different in the recursions for even and odd parity solutions. This 
implies that the cosine type series solution and sine series solution have different eigenvalues (Ae 2 ^ and Ao 2 ^ where 
m is the eigenvalue index) and they turn out to be complex conjugate to each other. The complex eigenvalue is a 
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consequence of non-Hermiticity of the Whittakcr Hill equation. These recursions can be written in a matrix form. 
The solutions can be expressed as an infinite trigonometric series and we have to truncate it at a point where the 
extra terms are not significant. We can obtain the eigenvalues and eigenvectors of the matrix which will provide the 
complete solutions to the Ince equations that are used in the solution to Whittaker Hill equations. Summarizing the 
angular solutions we can write 

oo 1 

Ae m {v) = e - lkRcoa »Y. a ™™s(r+-)v (52) 

r=0 

oo 1 

Ao m {v) = e~ ikR cos " Y, b mr sin(r + -)v (53) 

r=0 

Ae m (y) and Ao m (y) are solutions for the angular factor of the lower component of the separated spinor. Using the 
first order coupled equation (37 1 we can obtain the upper component angular eigenfunctions, 

Be m (v) = — - — {d v — ikR sin u)Ae m (v) (54) 

Bo m {v) = — (d v — ikRsmv)Ao m (v) (55) 

i\o m 

Be rn {u) and Bo m {y) are of the opposite parity to Ae m {y) and Ao m (v) due to the above operation and are given as 

1 °° 1 1 

Be m {v) = e - lfcflcos ^a mr (-r--)sin(r+-V (56) 

ZAC m Z Z 

r— (J 

oo 

Bo m (v) = — e^™ s "Vi mr (r+-)cos(r+-> (57) 

Note that for the sake of notation we always classify eigenfunctions according to the eigenvalues Ae^„ and Ao^. 
Functions corresponding to Ae^ get suffix "e" and corresponding to Xo 2 m get the suffix "o". 



B. Solutions to radial WHE 



Since we have obtained the eigenvalues by solving the angular equations, the eigenvalues can be used as parameters 
in the radial differential equations. 

The first method to evaluate the radial solutions is to replace i> —¥ i/i and (kR — > —kR) in the angular solutions 
(which is the same transformation that connects radial and angular WHE), the regular periodic (in solutions are 
denoted by Je and Jo. 



Je m (fi) = e kR cosh " Y c ™r sinh(r + (58) 

oo 

JoM = e lkR cosh " d mr cosh(r + -)n (59) 



2' 

r=0 



The lower component of the spinor can be obtained from coupled radial Eqs. (39) and (40). We denote the lower 
component radial solution with primes and keep this notation for all the lower component radial solutions. Note that 
prime does not imply derivative but is defined by the following operator acting on the upper component solutions. 



Je' m (fj) = —— {pix ikR sinh /it) J e m (^ti) (60) 
Jo' m {n) = -r—(d )M -ikRsixih.fi)Jo m (fi) (61) 

AO„ 
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The second linearly independent solution that is non-periodic in (i/x) is given as 

Ne m (fj,) = C: n {kR)^JeM + C c m {kR) fmr cosh(r + -) M (62) 

r=0 

oo ^ 

No m {n) = C° m {kR)^Jo m {ix) + C^{kR)Y,9 m r^Kr+ 1) )ii (63) 



r=0 

where the presence of factor of /i ensures the non-periodicity of the second solutions in i\x. Note that it is of opposite 
parity to the regular Je and Jo. This approach is similar to the calculation of non-periodic second solutions of the 
modified Mathieu equatioiM C^°(kR) are the normalization constants. 

The second method to calculate the radial solutions is the power series method. In this method we simply do a 
power series analysis in /i for the second order radial Ince equation ( |44| (note that such solutions do not capture the 
complex periodicity of the hyperbolic functions) . We then immediately get the two independent solutions of even and 
odd parity with the predetermined eigenvalues from the angular solutions acting as a parameter characterizing the 
different radial solutions. 

oo oo 

JeM = e lkR cosh " c r M 2r+1 , Ne m ([i) = e ikR cosh » £ f r y? r (64) 

r=0 r=0 
oo oo 

Jo m {n) = e ikR cosh " d r ^ N °M = z lkR c ° sh " E 9r M 2r+1 (65) 

r=0 r=0 

We are interested in studying the scattering problem which requires the radial solutions to have a well defined 
asymptotic form. We can evaluate the radial solutions as power series in ji and in series of sinh/u and cosh/z as 
described in the above mentioned methods. But these forms of solution diverge at large (x and therefore do not yield 
proper asymptotic forms. F ortun ately radial solutions to Whittaker Hill equations can be written as series of confluent 
hypergeometric functions^ * 20 * 23 ! which converge for all fi. We follow the procedure described in Ref. (jTSJ)- We start 
with radial Ince equation (44) 

a (/i) + 2ifci?sinh /ia (/x) — A 2 a(/i) = (66) 
and make the transformation z = cosh 2 ^. The resulting equation takes the form, 

z(z - l)a" (z) + (AikR z 2 - 4ikR z + z--)ol{z)- X 2 a{z) = (67) 

To extract the even and odd parity of solutions, we make the following functional transformations to the above 
equation. 

a(z) — \fz ae(z) (for even parity), a(z) = \J z — 1 ao{z) (for odd parity) (68) 
Making these transformations, we get the following differential equations for ae(z) and ao(z). 

z(z - l)ae" (z) + -(8ikR z(z - 1) + 4z - 3)ae (z) + (2ikR [z - 1) + - - Xe 2 )ae(z) = (69) 
z(z-l)ao"(z) + ~(8ikRz(z-l) + 4z - l)ao (z) + (2ikR (z) + ^ - \o 2 )ao{z) = (70) 

We have classified the eigenvalues as Ae 2 and Ao 2 for the even parity solutions and the odd parity solutions respectively 
(we know the eigenvalues from the angular eigenvalue Eqs. ( 48|49 1 and Eqs. ( 50|51 ). Solutions to the above equations 
can be expressed in terms of confluent hypergeometric functions 

oo 

ae m (z) = n M(n + 1/2, n + 2, -4ikR z) (71) 

n=0 
oo 

ao m (n) = E c ™ n M (n + 1/2, n + 2, -AikR z) (72) 

n=0 



10 



where the M are the regular hypergeometric functions satisfying the Kummer differential equation^ 

zM"{z) + (b- z)M'(z) - aM(z) = (73) 
The three term recursion relations for the coefficients c e m and c° m are 

c e m n ^AikR{n -\f + c e m „(n(n + 1) + 4ikR n + ~ - \e 2 m ) + (n + l)c e m n+1 = (74) 
c°m n -i^kR(n - i)(n + l -) + c° m n {n{n + 1) + MkR („+!) + !- Ao^) + (n + l)c° m n+1 = (75) 
And the full solution to the radial WHE (from the solution to the Ince equation) can be written as 



i — 1 

Jo m {p) = e lfcflcosh ^Jcosh 2 |^c^ n M(n+-,7i + 2,-4ifci?cosh 2 |) (76) 



2 

»=o 



JemOi) = e lkR cosh ^/cosh 2 5-l^ fl M(n + J, n + 2, -^fctfcosh 2 £) (77) 



2' ' 2* 

The second linearly independent solution can be obtained from the first solution using the following method. 

Fey m (n) = Je m (/i) J Je \^)2 d ^ ( 78 ) 

Mo 

Gey m {p) = Jo m {n) J Jq (^)2 rf ^ ( 79 ) 

Mo 

Fey and Gey are the second linearly independent solutions corresponding to Je and Jo. The lower component of the 



spinor can be evaluated by using the coupled equations ( 39 1 



1 

Ae r , 
1 

Ao„. 

(82) 



Fey' m (fj,) = ^^(8^- ikRsinhfi)Fey m (fi) (80) 
Gey' m (p) = (<9 M - ikRsinh /i)Gey m (p) (81) 



Solutions to the radial WHE for both upper component (ae m (/x), ao m (p)) and lower component (/3e m (/x), f3o m (fj,)) 
can be summarized in the combination of two linearly independent solutions as, 

ae m (n) = A e m Je m (n) + B e m Fey m (n) (83) 

/?e m (/i) = A e m Je' m (fi) + B^Fey'M) (84) 

ao m (n) = A° m Jo m (p) + B° m Gey m {p) (85) 

/3o ro (/i) = ^JoJJ/i) + 5^624^)) (86) 

A^j, and A° m , B^ are the undetermined coefficients and for notation sake we classify undetermined coefficients 
according to eigenvalues corresponding to Ae^ or Ao^. We identify undetermined coefficients with superscript "e" 
and "o" corresponding to the eigenvalues. A normalized choice for the undetermined coefficients would be 

Kn = A° m = zo^° m (87) 

BL = sin SL B°=sin8° 



Such a choice is helpful in formulating the scattering cross section in terms of phase shifts in the scattering amplitude 
with an being the phase shifts. 

Armed with all the solutions to the individual components of the separated Dirac spinor, we can now write the full 
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solution to the free Dirac equation as a superposition of all the eigenstates of the separated equations 

y,( M , v) = (ST) V ( r ^^^Kt^^ ( rl ) (89) 

m 



VI. EXPANSION OF INCOMING PLANE WAVE SPINOR IN TERMS OF WHITTAKER HILL 

EIGENFUNCTIONS 



The form of the incoming plane wave (see Appendix 



Eq. (Al8l) is given as e 



i k-'< 



COS : 

sin ! 



is the angle of 



incidence of the quasiparticle current). One of the requirements to construct the scattering cross section is to expand 
the incident plane wave in terms of eigenfunctions of the free Dirac equation which satisfy the following continuity 
condition (90 1 at [i = 0. 



(90) 



To write the plane wave expansion we take the following path. We write the free solution of the Dirac equation 
in elliptical coordinates as a linear combination of the eigenstates. One specific superposition of these eigenstates 
represents the plane wave spinor. Our aim in this section is to obtain these linear combination coefficients which 



represent the plane wave spinor. We can see that applying condition (90 1 to ip(/J,, v) is same as applying it to Y(fi, v) 
(see Eq. (32 )). This is because the transformation matrix which connects if) to Y cancels on both sides of the Eq. (p)0|. 



Y(Vi-v) Im=o= Y ifrv) 1/ 



(91) 



Applying the above condition, we can find the following constraint on the two undetermined constants (per eigenstate) 
appearing in the radial solutions due to the overall parity of the eigenfunctions. 



Bi 



B° = 



(92) 



Thus we see that the radial functions with constraint at ji = do not depend on Fey, Gey . Hence, the plane 
wave term only has Je and Jo terms which is analogous to the plane wave expansion in terms of Bessel functions in 
polar coordinates (see Ref. 28) which only contains regular J Bessel functions. Applying the appropriate boundary 
conditions on Y, we find the plane wave solution to be: 



(q T \Sr( (A e m Je m Be m (v) + A° m Jo m Bo m (v)) 
{ \ (i A e m J e' m Ae m (v) + i A° m Jo> m Ao m (u)) 



(93) 



The above solution with arbitrary constants A e m and A° m is an arbitrary superposition of eigenstates. We need to 
calculate the linear combination coefficients (as a function of 0) for which the expansion represents a plane wave spinor 
We write down the following expansion for the plane wave incident at an angle 8 with respect to the x-axis. 



gig i 



This step is important since we would like to control the angle of incidence of the incoming quasiparticle current. 



A k-r 



For details of the calculation of the plane wave expansion coefficients n e m and n° m refer to Appendix |D| 



Jo m Bo rn 



iJo m Ao m 



(94) 



VII. SCATTERING AMPLITUDE AND PHASE SHIFTS 



To study the scattering problem one requires well defined asymptotic forms of the full wavefunction spinor. The 
asymptotic form of the wavefunction determines the scattering amplitude in terms of the phase shifts of each partial 
wave with respect to the partial waves in the incident quasiparticle current. Thus, we consider the asymptotic limit 
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fi — > oo for the full wave function spinor ip = {ST )Y . 

K _ 1 ( cos | cosh f -sin f sinh %\sr^ ( ue m (n)Be m {v) + ao m (iJ,)Bo m (v) \ , Q _, 

W ' ' ~ (cosh n + cos t/) ^ sin f sinh f cos § cosh f J ^ ^ i(/3e m MAe m (i/) + /3o m ( M )v4o m (^)) J iyt,j 

To proceed further we need to use asymptotic forms for the radial functions. As /i — > oo, the asymptotic form (see 
Appendix [C|) for the radial functions is given as 

ae m (fi) ~ ^ -^e- 1 ^ {cos 5 e m + i sin 5 c m )e lh ^ e * (96) 

/?e m (^) - J : J-^e~ l i(cos5f n ~ i sin^Je"^ ^ (97) 

ao m (/*) - 7 ^^ e_if ( cos C + * sin eM (98) 

~ ^J^e-^Ccos^-isin^e-^ 6 " (99) 

For compact notation, we can write e M = ^ , where R is half the distance between the foci of the elliptical coordinates 
and r is the polar radial coordinate in the limit where the angular elliptical coordinate v approaches the polar angle 
<t>. In large fi limit elliptical coordinates reduce to polar coordinates. 



ae rn(r) ~ x /^- e -^e^ kr+s ^ (100) 



no.Jr) - x /-±- e ~^e« kr+S ^ (102) 



P ° m[r) ~ V^R e ' 4ie " (fer+C) (103) 

The asymptotic form of the ST transformation matrix can also be evaluated. As \i — > oo, cosh/i + cos v ~ cosh /i and 
cosh /i ~ sinh /j ~ ^- = . Applying these limits in the ST matrix of Eq. ( 30 1 gives 



s r ~ :r;^ i ci"-" 



We can then substitute the asymptotic forms of the radial solutions back into the full solution spinor of Eq. ( 89 1 to 
obtain 



7T7 *" 



(105) 

We would like to write the full wave function in the two suggestive parts requisite to set up the scattering problem 

ij) = l j J P lanewave _|_ ^scattered (106) 



The plane wave expanded in terms of the separated eigenstates Eq. (94 1 in the large /i limit takes the form 



r iane Wave = {ST) ^ J^ e -f £ e ikr (n e m Be m (9)Be m (v) + n° m Bo m {9)Bo m {v)) ( J ) (107) 



{ST) V^R e ~ li J2 e ~ ikri « Be ^ d ) Ae M +n° m Bo m (8)Ao m { V )) ( \ ) (108) 



13 



Now we construct the outgoing radial wave with appropriate asymptotic form. 



pattered = {ST) yj d ^ n ^ B e m (9)Be m (u)He m ( f i) + d m n° m Bo m (e)Bo m (v)Ho m (^)) ( J 



(109) 



He m (fj,) and Ho m (p) are the linear combinations of two linearly independent solutions to the radial WHE (Je m (/x), 



Fey m (fj,)) an d (Jo m (fi), Gey m (p)) which have the behavior of an outgoing radial wave (He m (r) ~ Ho m (r) 



) m 



the asymptotic limit. They play the role of Hankel functions^ (linear combination of Bessel J and Bessel Y) which 
appear in the study of scattering problems in polar coordinates. d e m and d° m are the undetermined coefficients. The 
asymptotic form of the scattered wave is given as 



scattered 



(n -> oo) = (ST) 



irkR 



< > (d e m n e m Be m (6)Be m (v) + d m n° m Bo m (6)Bo m (v)) 



(110) 



Now we compare our full wave function spinor, Eq. (89 1 and with the wave functions written in a suggestive form 
in the asymptotic limit. e lkr and e~ lkr multiplied by the angular functions in v are independent functions. Hence 



corresponding to the four independent angular functions. 



we can equate their coefficients in Eq. (1061 and Eq. (105 1. We get four equations for the undetermined coefficients 



= (l + d e m )n e m Be m (9) 

= (l + d° m )n° m Bo m (9) 

= (n e m )Be m (6) 

= (n° m )Bo m (d) 



(111) 
(112) 
(113) 
(114) 



Solving the above four equations we can write d e m and d° m (which are the undetermined coefficients of the scattered 
wave) in terms of the phase shifts 6^ and 6° t . 



d e m = (e 2lS - ~l),d° m = (e 



2i&° 



1) 



(115) 



To write the full form of the scattering amplitude we multiply the scattered wave by the ST transformation matrix 
in its asymptotic form. 



1> 



scattered 



J2(d^n^Be m (e)Be m (u)+d m n° m Bo m (e)Bo m (p))^ < ( '" s * 



(116) 



cos 
sin - 



represents the quasiparticle current going in the radial direction (see Eq. ( A21 1). The scattering amplitude 



can be extracted from the asymptotic form of the scattered wave and is given as 

T 



l)n e m Be m (6)Be m (v) + (e 



1) n° m Bo m (6)Bo m (v)) 



(117) 



Hence, we were successful in constructing the scattering amplitude (analogous to the general form of scattering 
amplitude (Eq. ( |A20 )) in elliptical coordinates. The only thing that remains is to calculate the phase shifts ^ and 
5!^, and for that we need to impose the conditions for the branch cut on the full wave function spinor, Eq. (89 1. 



VIII. SCATTERING CROSS SECTION WITHOUT BRANCH CUT (BERRY PHASE PARAMETER 

B=0) 



Before going on to the case with the branch cut, we make a quick check on our scattering amplitude for the case of 
no branch cut or B=0 in Eq. 0. We expect this trivial case to yield no scattering of quasiparticles. For the case of 
no branch cut between the foci of the ellipse (/i — 0), we impose the following condition on the wave function spinor 



(118) 
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At n — 0, we have the following behavior for Je m (fi), Je' m (/i), Jo m (/x) and Jo' m (fj,) 



Jo m (0) ^0, Je m (0) = (119) 
Je^(0)^0, J O ' ro (0) = (120) 



At fi = 0, the second independent solutions are all nonzero, 



Fey m (0) £0, Gey m (0) + (121) 
Fe y ; n (0)^0, Gey' m (0)^0 (122) 



Applying condition (1181, all the terms containing even parity angular eigenfunctions cancel out and we can write the 
remaining terms as 

(cos S e m Je m (0)+ sin 5 e m Fey m {Q))Be m {v) \ f 
i (cosd° n Jo' m (0) + sxn8° m Gey' m {0))Ao m {v) ) \ 

Substituting for the values of the radial functions at [i = 0, we have 

(sin5^Fey m (0))Be m (v)) \ ( 
i (s\n5° m Gey' m (0))Ao m {v) ) \Q 

From the boundary condition at the origin we get the following constraints on the undetermined coefficients, 

sin<5^ = 0, sin<5^ = (123) 

The above expression gives phase shifts as 6^ = and 5° n — 0. Putting the obtained phase shifts back into Eq. ( |117 l 
for the scattering amplitude, we obtain 

f(e,v)=0 (124) 

Hence, we recover our trivial result that without the branch cut (and without superflow) there is no scattering. Now 
we move to the interesting case of quasiparticle scattering with the branch cut. 



IX. SCATTERING CROSS SECTION DUE TO A BRANCH CUT (B=l) 

For the case with a branch cut between the foci of the ellipse (/i = 0), we set the Berry phase parameter B=l which 



sets the condition imposed on the wave function spinor Eq. ( 89 1 to 



V>G« 5 -v) U=o= -ipifi, v) | M=0 (125) 

Applying the above condition and using the values of the radial eigenfunctions at \i = 0, all the terms containing odd 
angular eigenfunctions cancel out and we can write the remaining terms as 

(cos 5° m J o m (0) + sin 8° n Gey m {0))Bo m {v) \ f 
i (cos^Je^(O) + smS? n Fey' m (0))Ae m (v) ) \ 

From the boundary condition at the origin we get the following constraints on the undetermined coefficients, 

cos^Jo ro (0) = -sin^Gey Tn (0), cos^Je^(O) = -mi6^Fey' m (0) (126) 

In other words 



Fey' m (o) 

Jo m (0) 
Gey m (0) 



tan<5°„ = - i "'^, (128) 
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Since tan 5^ and tan 5^ are complex conjugates of each other we have to account for the relative sign between 8^ 
and <$£i while calculating the inverse tangent in the above relation. Substituting for JfL and 5° in the definitions of 



the phase shifts, Eq. (1151 we obtain 



(e 2 ^- - 1) 
(e 215 - - 1) 



-2Je' m (0) 



Je' m (0)+iFey' m (0) 

-2Jo m (0) 
Jo m (0) - i Gey m (0) 



(129) 
(130) 



Using Eq. ( 129 1 and Eq. ( 130 1 we can completely evaluate the scattering amplitude and the differential cross section 

Jo m (0) 



for quasiparticle scattering due to branch cut without the superflow 

Je' m (0) 



Yd 6 U £ ( Je' 



-)n e m Be m {9)Be m {v) + (- 



,(0)+!Fe<(0)' , ' m "'"" v "'"'"" v "' ' v Jo m (0) -*Gey m (0) 
We can write the exact differential cross section for the quasiparticle scattering in terms of cp = v — 9, 



)n° m Bo m {9)Bo m {v) (131) 



da _ 4 
dp irk 



Je' m (Q)+iFey' m (0) 



)n e m Be m {9)Be m (p + 0) 



Jo m (0) 



Jo m (0)-i Gey m (0) 



)n° m Bo m (0)Bo m (p + 0) 



(132) 



Recall that at the outset of this calculation, we shifted the origin of momentum space to the center of node 1. Thus, 
in the discussions that followed, we have been considering quasiparticles scattered from one state in the vicinity of 
node 1 to another state in the vicinity of node 1. The resulting cross section is therefore only the cross section for 
these node-1 quasiparticles. However, given a quasiparticle current in any particular direction, quasiparticles from all 
four nodes will contribute equally. Thus to obtain the physical cross section, we must average over the cross sections 
for quasiparticles at each of the four nodes. Our results for node 1 can be easily generalized to node j = {1, 2, 3, 4} by 
transforming coordinates to those appropriate to node j. In accordance with the d-wave structure of the gap, we can 
define a local coordinate system at each of the four nodes with a ki axis pointing along the direction of increasing 
efe and a k2 axis pointing along the direction of increasing Note that while nodes 1 and 3 define right-handed 
coordinate systems, nodes 2 and 4 define left-handed coordinate systems. We can therefore transform from node 1 to 
node j simply by rotating our incident and scattered angles (0 and v) and then changing the sign of these angles to 
account for the handedness of the local coordinate system. 



Node 1: 


1 = 


Node 2: 


02 = -{9- 


Node 3: 


9 3 = e + TT 


Node 4: 


04 = -{9 + 



I/3 = V + 7T 



<P1 


= v\- 


Hi 


= V 




= v 2 - 


02 


= -V 






03 


= f 


V 3 -! 


= 1*4- 


04 


= -V 



(133) 



Thus, to obtain results for quasiparticles about node j, we need only input each 9j and take the output as a function 
of ( — iy +1 Lp. Then the physical cross sections are 



do 

dp 



4 ^ 



da 

dp 



(134) 



X. DIFFERENTIAL CROSS SECTION RESULTS FOR BERRY PHASE SCATTERING 



In this section we plot the differential cross section for quasiparticle scattering from branch cut (no superflow) for 
several cases. We define the distance between the two vortices by the dimensionless parameter kR. Incident angle for 
the current is described by as shown in Fig. [2] The thick dots on the foci of the ellipse depicts the vortex cores and 
the thick line joining the cores denotes the branch cut and is also the fi=0 line in elliptical coordinates. The wiggly 
arrows represent the incident quasiparticle current. The angle of incidence is with respect to the x-axis. The plane 
wave spinor representing quasiparticle current is incident on the vortex cores and part of it acquires a Berry phase 



1G 



factor of (-1) between the vortex cores. The scattering contribution is entirely due to this effect. 



v=jc/2 




v=0 



Figure 2: The above picture depicts the scattering of quasiparticles due to the Berry phase effect. The Berry phase effect is 
denoted by the finite branch cut shown by the thick line joining the dots. The dots represent the vortex cores coinciding with 
the foci. Wiggly lines denote the incident quasiparticle current. 6 is the incident angle of the quasiparticle current with respect 
to the x-axis. 



kR=0.1 













6/f v jrk 7 





dtr / 4 \ 
die l irk ' 




dcr { _4_\ 

iif ^ jrk 





Figure 3: In the above figure vortex cores (foci of ellipse) are depicted by dots and the line joining them is the branch cut. 
Vortex cores are separated by dimensionless length kR=0.1. 8 is the angle of incidence of the quasiparticle current. For small 
inter-vortex separation, the ellipse looks like a circle which indicates near-circular symmetry in the scatterer. We plot single 
node differential scattering cross section for quasiparticle current incident at different angles 6. The plots of the scattering cross 
section emphasize the near-circular symmetry with respect to the incident angle 6 due to small inter vortex separation. 



For very small inter vortex separation (kR=0.1) the ellipse looks like a circle (see Fig. |3J) . In this limit we expect 
near-circular symmetry in the plots of differential cross section. We see that the quasiparticles see relatively small 
branch cut which results in the differential cross section being almost independent of the incident angle 6 (see Fig.|3|. 
Note that this is not a very good limit physically since we can no longer ignore the presence of other vortices in the 
sample. When we stretch the vortex cores apart (kR=l) the scatterer becomes more elliptical (see Fig. |4|. This is 
reflected in the elliptical symmetry we see in the cross section plots as we rotate the incident angle of the quasiparticle 
current (see Fig. [4]). We see that we get the same plots for differential cross section if we rotate the incident angle by 
7r, which reflects the symmetry of the scatterer (symmetric under 7r rotation) . 
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Figure 4: In the above figure vortex cores are further apart with dimensionless length kR=1.0. With the increase in inter- 
vortex separation the scatterer becomes more elliptical and plots show expected elliptical symmetry in the single node differential 
scattering cross section. Also note the increased magnitude of scattering cross section which can be attributed to the increase 
in the length of branch cut. In other words, more quasiparticles hit the branch cut 




Figure 5: In the above figure vortices are further apart with kR=3.0. The plots of the single node scattering cross section show 
elliptical symmetry. We also see increase in the magnitude of scattering cross section as compared to the case of kR=1.0 

As we increase the inter-vortex separation further to kR=3.0, we obtain the case of a highly elliptical scatterer. For 
this case the magnitude of the scattering cross section increases as compared to the case of kR=l (see Fig. [5]). We 
observe the scattered current sweeping closer to the forward direction for higher kR. For kR=l and kR=3, we see 
maximum scattering for the case of 6 — n/2. At this angle the quasiparticle current is normally incident on the branch 
cut and results in maximum exposure to the Berry phase effect. Mirror symmetry about 9 = ir/2 is seen in the cross 
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section plots. Thus our results for the scattering cross section are consistent with the geometry of the scatterer. We 
note that the scattering is reflectionless, or in other words, there is no backscattering of the quasiparticle current. 
Absence of backscattering due to the Berry phase has been previously reported in the literature for the case of carbon 
nanotubes (see Ref. [55). 

We must now average over the scattering contribution due to quasiparticles from all four nodes. The four node 
average has been performed as prescribed in Eq. (1331 and Eq. (134). After averaging over four nodes, we still see 9 
dependence in the differential cross section (see Fig. |6 I . The resulting cross section is n/2 periodic with respect to 9. 
This is a consequence of the definition of 9 at each node (see Eq. (133)). We should keep in mind that in the setup 
we consider (cuprate sample), the pair of vortices are not always aligned along the x-axis as shown in Fig. [2] 




kR=J — FOURNODE AVERAGE 





4; ' .Tfc 1 






Figure 6: Four Node Average Differential Cross Section. We plot differential cross section averaged over the contributions 
of quasiparticles from all four gap nodes, for quasiparticle current incident at various angles 6 and for inter-vortex separation 
kR=3.0. Results are 7r/2 periodic with respect to 9. 




:.:R 



Figure 7: ka\\ plotted versus increasing inter vortex separation, averaged over the incident angle 9. The solid curve shows the 
transport cross section for the Berry Phase scattering case. Dashed curve shows the transport cross section for the superflow 
scattering. Inset shows kau plot for the case of superflow scattering of quasiparticles plotted for very high kR values^. 

In Fig. [7] we plot the total transport cross section, fc<7|| as a function of inter-vortex separation kR. 

dip — (1 — cos <p) (135) 
We notice that the transport cross section goes to as kR— ^0. This is expected since the Berry phase effect (branch 
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cut) is negligible for very small values of kR. With the increase in kR, the longitudinal cross section increases rapidly 
and then saturates for kR>l. On the basis of the transport cross section plots obtained for the case of superflow 
and Berry phase scattering, one can make an intuitive comparison of these two effects. In the superflow paper^l, we 
have neglected the Berry phase effect by applying periodic boundary conditions to the quasiparticle wave functions. 
This amounts to scattering of quasiparticles from vortices with superflow potential with the strength of two vortices. 
To calculate the transport cross section for the case of a branch cut between the two vortices, we have neglected the 
superflow contribution. Hence we have neatly isolated the scattering contributions due to these two effects, which 
gives us an opportunity to compare these two effects. Before such a comparison, we must treat the case of the Berry 
phase scattering on an equal footing with the superflow scattering. Due to the two-center nature of the regularized 
Berry phase effect, we are dealing with elliptical geometry in this case. Upon performing the four-node average, we 
see that this elliptical symmetry shown in the differential cross section plots (See Fig. [5| has been reduced to near 
circular symmetry even for the highly elliptical case of kR=3 (Fig. [6| . Also, since there is no preferred orientation 
of the branch cut, one can average over the alignment of the branch cut with respect to the x-axis. The final DCS 
averaged over this alignment will have no elliptical symmetry or skew scattering. Based on the above arguments, we 
may directly compare the transport cross section due to superflow potential of two vortices of radii kR on top of each 
other to the Berry phase scattering due to two vortices separated by distance kR (averaged over incident angle 9). In 
both cases, kR parameterizes the dimensionless energy of incident quasiparticles. kR also determines the size of the 
vortex for the case of superflow scattering and the length of the branch cut for the case of Berry phase scattering. For 
the superflow case, we see a steep increase in the transport cross section (after averaging over 4 nodes) followed by 
saturation for kR > 5. For the case of Berry phase scattering, the increase in transport cross section is steeper than 
for the superflow case. For the Berry phase scattering the magnitude of kern saturates for kR > 1. This shows that 
the Berry phase is the more important effect of the two for kR < 5 (high field low temperature regime). Magnitudes 
of transport cross sections for higher kR (weak field high temperature regime) are similar for both the superflow and 
Berry phase processes. On the basis of the plots of transport cross section for both cases, one can conclude that 
the transport cross section due to the branch cut dominates for kR < 5, and is of similar order to the superflow 
contribution for kR > 5. 



XI. CONCLUSIONS 



In this work, we calculated the Berry phase contribution to the scattering of quasiparticles from vortices in a d-wave 
superconductor. We simplified the Bogoliubov-de Gennes Hamiltonian by applying a singular gauge transformation. 
This transformation extracts the phase from the gap function and encodes it in the antiperiodic boundary conditions 
imposed on the wave function. Within the single vortex approximation, this antiperiodic boundary condition (Berry 
phase) manifests as a semi-infinite branch cut such that with each trip around the origin, the wave function changes 
sign. We neglected the superflow contribution and considered the scattering of quasiparticles due only to the presence 
of this antiperiodic boundary condition. We found the scattering cross section for this case to be divergent in the 
forward direction. In order to regularize the Berry phase effect we considered the two vortex problem (two Aharonov- 
Bohm half fluxes) which has a finite branch cut between its cores. To solve this two-center problem, we chose to 
work in elliptical coordinates. Elliptical coordinates provide an advantage in implementing the branch cut condition 
on the wave function spinor. We can turn on the Berry phase effect by simply imposing the boundary condition 
(see Eq. Q) on the wave function spinor. We separated the (2+l)d Dirac equation in elliptical coordinates and 
found that the separated equations were Whittaker Hill equations (WHE). We solved the eigenvalue problem for 
the angular WHE. With the calculated eigenvalue as a parameter, we obtained the two linearly independent radial 
solutions with well defined asymptotic behavior. We developed a plane wave expansion for the incident quasiparticle 
current in terms of the separated solutions of the WHE. Using a partial wave analysis, we expressed the scattering 
amplitude and differential cross section in terms of phase shifts. We obtained the phase shifts by imposing the branch 
cut condition (Berry phase effect) on the full wave function spinor. We analyzed the scattering cross section due to 
the Berry phase effect for different separations of the vortex cores. We have also presented the variation of the total 
transport cross section as a function of the inter-vortex separation. We have also given a qualitative comparison of 
transport cross sections for the Berry phase and the superflow scattering mechanisms in Sec.|X] Berry phase scattering 
of quasiparticles discussed here is not restricted to the case of d-wave superconductors. With some modifications of 
the incident plane wave, our problem becomes that of general relativistic scattering in two dimensions due to two 
Aharonov-Bohm half fluxes. In this work, we have neglected the superflow contribution. Single vortex scattering 
due to a circulating superflow is considered in a separate papei^. The problem that we have considered in both the 
manuscripts is a simplified version of a more complicated scattering process. We have made a series of approximations- 
to tackle the problem in its simplest form. This work should be treated as a first step forward to understand the 
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complicated and important issue of scattering of quasiparticles from vortices. Deviating from the analytical setup, 
we can improve the model by considering the anisotropic Dirac spectrum, internodal scattering, and an even more 
rigorous description of plane waves. All these effects become important once we move away from the weak field limit . 
The linearized version of the BdG equation is limited to low energy quasiparticle excitations and one must use the 
full BdG equation and solve it numerically for higher energy cases. One way to include the above mentioned effects 
is to consider a sea of vortices in a lattice model. Such a vortex lattice calculation has been considered by Melikyan 
and TesanovicM In the appendix of their paper, they have set up the two vortex scattering problem in elliptical 
coordinates. In this work, we have successfully obtained the exact scattering solutions to the two vortex problem 
in elliptical coordinates. The results of our calculations along with the vortex lattice calculations provide a greater 
insight into the bigger picture of quasiparticle scattering from vortices. To this end, our simplified model of double 
vortex scattering in an analytical framework is an important result. Our next step will be to consider both effects 
within the double vortex model by including the superflow that circulates around the vortices in the presence of the 
branch cut that lies between them. We expect that this analysis, left for future research, will provide insight about 
not only the relative importance of the two contributions but also the interference between them. 
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Appendix A: Bogoliubov-de Gennes Equation 

The setup for this problem is described in detail in Ref. [5] We summarize here. Consider the Bogoliubov-de Gennes 
(BdG) equation for a d-wave superconductor in the presence of a constant perpendicular magnetic field, A = ^Hrcj), 
and with an order parameter that winds once about the origin, A(r) = A e l( ^: 

H'* = E* H '=(§* ( A1 ) 
H' e = ^( P --a) 2 -E F (A2) 



2m V c 

A' = ^{p x , {p y , A(r)}} - A- A(r)(W) (A3) 

Here p = —ihV, {a, b} = (ab + ba)/2, and E is the quasiparticle energy. The form of the gap operator enforces 
the ci-wave symmetr y^ * 31 !. Upon circling an hc/2e vortex, the quasiparticle acquires a Berry phase factor of (-1). 
This fact is encoded within the complex differential form of the gap operator, A. We simplify the Hamiltonian by 
effectively stripping the gap function, A(r), of its phase. This is done by applying the singular gauge transformation 

U = J /2 *(r) = U-H(t) H = U~ l H'U, (A4) 



= E<£ (A5) 



e 1 ^ 2 

In this gauge, known as the Anderson gauge, 



H = r 3 ^[(p + r 3 P s ) 2 - p%] +n^- [2 PxPy ] (A6) 

IpF ^PF 
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where 

p. w _»v#-Ja -!(!-£)# ( A7) 

is the gauge invariant superfluid momentum (superflow), vt = pf/tu, V2 = Aq/pf, and R = y/hc/eH . In effect, 
the Berry phase contribution has been extracted from the Hamiltonian and encoded in the antiperiodic boundary 
conditions imposed on the wave function. While the original wave function was defined with periodic boundary 
conditions, ^(r, 4>) = ^(r, <f> + 2tt), the transformed wave function is not single-valued and has antiperiodic boundary 
conditions, <&(r, <j>) — — $(r, cj) + 2ir). Hence, we have introduced a branch cut such that with each trip around the 
origin, the wave function changes sign. Note that the Berry phase effect is not a consequence of the choice of gauge. 
We have used the singular gauge transformation to extract the Berry phase contribution from the Hamiltonian and 
encode it in the boundary conditions of the wave functions. By definition, all observables, such as differential cross 
section or transport coefficients, are independent of this gauge choice. For other problems, other gauge choices are 
optimal. A nice discussion of this is provided for the case of vortex lattice by Franz and TesanovicP and Vafek et 
aP. We can further simplify our Hamiltonian by shifting the origin of momentum space to the location of one of the 
nodes. Shifting to node 1 



PF +Px Py Py (A8) 



we find that 



H = H D + H C (A9) 
H D = v f [p x r 3 + ap y Ti + P sx ] (A10) 

H C = ^- [ip 2 + P?)t 3 + 2P S • p + a2p xPy n] (All) 

ZpF 

where a = Vijvj and we have used the fact that P s — P s (r)(j> to commute p with P s . Here H is written as the 
sum of a linear (Dirac) Hamiltonian, Hjj, and a quadratic (curvature) Hamiltonian, He- The second (curvature) 
term is smaller than the first by a factor of E/Ep- We will focus on the dominant term, Hp,- In order to study the 
quasiparticle scattering from vortices, we must consider the nature of quasiparticle current in a d-wave superconductor. 
Since the incident and scattered currents will be considered in the far field where the quasiparticles are free, we wish 



to determine the quasiparticle current as a functional of $ for P s = 0. Setting P s = in Eq. (A6l we find that the 
BdG Hamiltonian becomes 

H=( I ! e H e = -^V 2 ~E F A = -^2d x d v (A12) 

V A* -H* J 2p F 2p F x v y ' 

Following Re fs. [5] and [331 we can write down the quasiparticle current corresponding to the BdG Hamiltonian. Once 
again, it is convenient to shift the origin of momentum space to a nodal point. Shifting to node 1 yields 

j = 3d+3c (A13) 
3d = v f & fait + any)® (A14) 



V f T 

Jc = — Im 



+ d® i d® 

$ T (r 3 x + ariy)— + <£> T (T 3 y + an*) — 



(A15) 



PF 

where a = v^/vf. In what follows, we will focus on the dominant term corresponding to the Hamiltonian H/j. 

To proceed, we must obtain a general form for the scattering cross section. We consider a plane wave, with quasiparticle 
current in the incident direction, scattering off a vortex as a radial wave, with quasiparticle current in the scattered 
direction. If the incident momentum is k = [k, 9) and the final momentum is k' = (fc, <p), then the incident direction 
is the direction of the group velocity at momentum k and the scattered direction is the direction of the group velocity 
at momentum k'. For general, anisotropic nodes, the group velocity need not be parallel to the momentum. However, 
for the isotropic case that we consider 

v G (k) = = v f ^-x + v 2 ^ L y = v f (cos 0x + sin 9y) = u/k (A16) 

OK -C/fe -C/fc 



and the group velocity and momentum are parallel. Therefore, if <F; denotes the incident wave function and $ s denotes 
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the scattered wave function, then we require 

Jd ~ (cos flx + sin 9y) ~ k jo [$ s ] ~ (cos 0x + sin </>y) ~ k' ~ f . 



(A17) 



Inspection of the form of the current functional, jn = Vf&' (T3X + riy)$ reveals that the appropriate incident plane 
wave is 



COS ; 

sin ! 



(A18) 



Note that outside the vortex, quasiparticles are subject to neither an order parameter phase gradient nor a magnetic 
field. Thus, the incident wave function is a plane wave. This is consistent with the well-known results of Franz and 
TesanovicP who showed that the low-energy quasiparticle states of a d-wave superconductor in the vortex state are 
Bloch waves of massless Dirac fermions rather than Landau Levels. (For a discussion of the analysis that let to this 
important result, the reader is referred to Refs. 7 34 36l) 



cos 2 - — sin 2 - I x 
2 2 . 



2 sin - cos - I y 
2 2 ' 



u/k. 



(A19) 



Note also that this form solves the BdG equation, as it must in the absence of the vortex. The appropriate scattered 
radial wave is then given by 



*.(r)= 6*5/(0,0) 



e tkr 1 cos- 



y/r \ sin 



(A20) 



Iff 



cos" sin" — I x + I 2 sin — cos — I y 

2 2 V 2 2 ' " 



v f — r - 



(A21) 



Here f(</>, 9) is the scattering amplitude and the e 1 ^!" 1 prefactor has been added to make the wave function single- 
valued. 



Appendix B: Separation of Variables for massless spin-1/2 2D Dirac equation in elliptical coordinates 



We present here the details of the steps (see Sec. IV I leading to the separation of variables in elliptical coordinates 
of a two-dimensional Dirac equation for massless spin-1/2 fermions. This calculation is in the spirit of work done in 
Ref. [T5] An alternate method can also be found in Ref. [T3J The 2-D Dirac equation is given by 

[ 1 O d t + 1 1 d X + 1 2 dy}^ = (Bl) 

where the 7's satisfy the following anticommutation relations 

{ 7 Q , 7 ^} = 2 5 ^ (B2) 

with Minkowski metric given by g = diag(l,—l, —1). In order to rewrite the Dirac equation in curved coordinates, 
let us introduce the following coordinate transformation 

x = .f{^,v), y = g{n,v), t = t. (B3) 

The reader should note that the transformations are kept general and the choice of elliptical coordinates will be made 
when required. The general transformations given by /(/i, v) and g(p, v) must satisfy the condition that / + ig is 
holomorphic or complex differentiable in the u + iv plane which leads to the following Cauchy-Riemann equations, 

U =9v, 9 l i = -fu (B4) 



23 



Using Eq. ( B3 1 and Eq. (B4| in Eq. (Bl). one can easily write the dirac equation in curved coordinates as 



where 



[fd t + ^ + ^—d v ]^ = 



7 1 = ^if^-Ul 2 ) 

7 2 = ^W + /„7 2 ) 



with the Lame Metric given by h = < ft + g*. We introduce the following transformation matrices: 



1 , * 



I A i I 1 2, 



with 4> — arctan(j a ). We use Eq. (B8) in transforming the Dirac equation. That is we perform 



which results in 



1 2 

[d t + 1 -^ L d^ + 1 —^d v ]^ = Q 
h h 

where the transformed spinor satisfies 

$ = s _ V 

To separate the time variable we introduce the following operator definitions. 

k 2 = d t , ki = (^d„+^d v ) 



The Hamiltonian Eq. (BIO I will read 



[ka + h}$ = 

with fci and hi satisfying the commutation relation 

\ki$x\ = 

To separate the time variable, we introduce k such that 
which immediately gives 

Here, we make the following choice for a two-dimensional representation of the Dirac matrices, 

7° = T 2, 7 1 = * T i i 7 2 = -* T 3 



In this representation, Eq. (B16) reads 



[nd ll + Tid v -ikh]Q = Q 



(B5) 

(B6) 
(B7) 

(B8) 

(B9) 

(BIO) 
(Bll) 

(B12) 

(B13) 
(B14) 
(B15) 
(B16) 
(B17) 
(B18) 
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The presence of h in Eq. ( B18 ) forbids us to write Eq. ( B18 1 as sum of two commuting differential operators. Therefore 

(B19) 



we will introduce a similarity transformation T(/i, v) acting on the Dirac operator and the spinor 

T = e e iaT2 
with 



To be more explicit, we do the T-transformation on Eq. ( B18 1 in the following way 



which after some algebra gives, 



with 



T[T 3 d„ + nd v - ikhlTT- 1 ® = 
[r 3 a M + T X d v - ikhe t2aT2 ]Y = 



Y = T- 1 ® 



It should be noted that because a([x, v) in Eq. (B19) is arbitrary we choose it to be of the following structure, 

i2ar 2 _ [Op) + ih { v ) T 2\ 



(B20) 
(B21) 
(B22) 
(B23) 

(B24) 



This specific form of the transformation matrix T cancels the factor h (Lame metric) in the Dirac equation which 
mixes the /i and v variables. It is trivial to check that a(/i) 2 + b(v) 2 — h 2 and the structure for /3(ji, v) can be obtained 
using Eq. ( B20 ) . At this point we move to elliptical coordinates with Eq. ( B3 1 taking the form 



From (B25 B26l we get 



x = /(/x, v) = R cosh fi cos v 
y = g([A, v) = R sinh/i smu 



h 2 = R 2 sinh 2 \i + R 2 sin 2 v 



Using above (B27l along with the fact that a(/x) 2 + b{v) 2 = h 2 gives us by comparison the following 

o(/Lt) = R sinh fi, b{v) — R sin v 
Now we introduce the following operators 

L\ = r 3 9 M - ika(fi) 
L2 = T Z d u - ikb{y) 

y = i^+i^r^z 

Using Eqs. (B29l, (B30), and ( B31| ) we can finally express the Dirac equation Eq. (B22l as 



P 2 - ikr^a^ + k 2 a 2 ) + (<9 2 - ikr 3 b u + k 2 b 2 )]Z = 



Defining 



(B25) 
(B26) 



(B27) 

(B28) 

(B29) 
(B30) 
(B31) 

(B32) 
(B33) 
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and introducing a separation constant A gives us four 2nd order ordinary differential equations. 

A 2 )«( M ) 

fcV - X 2 )/3(p) = 



(d 2 - ika^ + fcV - A 2 )a(» = 



(d 2 + ika„ 



) 2 
ii 

Ol 

Ol 



ikby 
ikb,. 



■k 2 b 2 
k 2 b 2 



X 2 )A(u) 
X 2 )BM 



After putting back a — i?sinh/i, a M = Rcoshfi, b = Rsvctv, and b v = Rcosv we obtain 

{d 2 u -ikR cos v + k 2 R 2 sai 2 v + \ 2 )A{v) = 
(dl +ikRcosv + k 2 R 2 sin 2 v + \ 2 )B(v) = 



(9„ — ikRcosh/i + k 2 R 2 sinh fi — A 2 )a(/Lt) 







(<9 2 +ifci?cosh^ + fc 2 i? 2 sinhV- A 2 )/%) = 



(B34) 
(B35) 
(B36) 
(B37) 



(B38) 
(B39) 
(B40) 
(B41) 



which are equivalent to 4 coupled first order equations that connect the upper and lower components of the wave 
function spinor 



(d v - ikRsinv)A(v) = i\B{v) 

(d v + ikRs\nv)B(v) = i\A(v) 

(dfj, — ikRsinh fi)a(fi) — A/3(/x) 

(d^ + ikRsmhfi)/3(fj,) = Xa(^i) 



(B42) 
(B43) 
(B44) 
(B45) 



Thus we have reduced the 2D massless Dirac equation in elliptical coordinates to a problem of four decoupled ordinary 
differential equations. These separated radial and angular equations are known as the Whittaker Hill equations 
(WHE). 



The explicit form of the transformation matrices can be evaluated. Using Eq. (|B31|) we have 

Y 



iP(p)A(u) 



We write below a more transparent and dimensionless form for S (Eq. (B8l) and T (Eq. (B19l) 

S = 



'R 
T = 



cos ■ 
sin ■' 



sin ■ 



cos 



2 2 

cos a sin a 
- sin a cos a 



which immediately gives ST as 



st = Vr 



cos(a - <j)/2) sin(a - (f>/2) \ 
— sin(a — <fi/2) cos(a — <j>/2) J 



In elliptical coordinates, a, /3, and 4> take the following form, 



1 _-, . sin v 
a — — tan ( -— — ) 



v sinh /x' 



1 



. _1 , cos V . 

= --tanh ( — ) 

2 cosh /i 

tan _1 (coth \i tan v) 



(B46) 

(B47) 
(B48) 

(B49) 

(B50) 

(B51) 

(B52) 
(B53) 
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and we can rewrite j3 using the log formula for tanh 1 

P = ln( 



COSh [L + COS ^_i/4 



cosh /i — cos v 
eP_ _ 1 

\J (cosh ^ + cos v) 



After using some trigonometric identities, these expressions for a and satisfy the following relations 

4> tan_1 (sar7l) ~ tan-^coth yt tan i/) 
""2 = 2 

2(o-^) = tan-H - sil ^f^ ) 
2 cos i/ cosh [i + l 



(B54) 

(B55) 
(B56) 

(B57) 
(B58) 



To evaluate ST we need to calculate its components in terms of elliptical coordinate variables. Tangent, sine and 
cosine of the double angle can be evaluated as 

— sin v sinh [i 



tan(2(a — 


P 


sin(2(a — 


2 n 


cos(2(a — 


2" 



COS V 


cosh [i 


+ 1 


— sin 


v sinh 


[i 


cos v 


+ cosh 


f 


cos v 


cosh [i 


+ 1 



cos v + cosh [i 



Using the half angle formulae we can finally evaluate the components of the transformation matrix in terms of the 
coordinate variables. While reducing the components to half angle there is arbitrariness in the choice of the sign of the 
trigonometric functions in each quadrant. We make a choice of signs such that our final transformation matrix is 4-7T 
periodic in angular coordinate v. We could as well make a choice which would result in 2~k periodicity of transformation 
matrix. Any choice can be compensated by choosing the periodicity of the angular solutions accordingly. For example, 
we choose 4-7T periodic eigenstates for the angular WHE (see Sec. |v| . 



sin(a 1 

v 2' 



cost a 

v 2 1 



sin | sinh ^ 
-^/cos v + cosh [i 
cos | cosh t£ 
Vcos v + cosh Jt 



(B59) 

(B60) 
(B61) 



Using Eqs. ( |B49J flB55| , ( |B59| , and (JB60J) along with the fact that ip = ST Y, we have 



1 



(cosh [i + cos v) 



cos | cosh ^ 
sin | sinh ^ 



— sin | sinh ^ 
cos | cosh ^ 



a{[i)B{v) 
if3([i)A(v) 



Appendix C: Asymptotic form for radial solutions 



To study the scattering cross section, we have to calculate the phase shifts of the scattered wave function at large 
radial distances where the detector is placed. Hence, the asymptotic form of the radial solutions is an important piece 
of information in setting up the scattering cross section. In this section we give the asymptotic form of the first and 
second radial solutions to WHE. In the Sec. |VB| we expand radial solutions in terms of confluent hypergeometric 
functions (CHF) (see Eqs. (69 1 and ( 70 1 ) . We now write the asymptotic form of these solutions. As /i — > oo in the 



27 



series of CHF, only the leading term in the sum contributes 



Jo m (/i) « jkRcoshu CQsh | M (I ; 2 , -4ifci? cosh 2 | ) (CI) 
Je ro (/i) « e lfcflcosh ^ sin h| M(i,2,-4ifci?cosh 2 |) (C2) 



The asymptotic form of the CHF is well knowrPS and is given as 



M(i 2,-4i/ci?cosh 2 « — L./ — =— (C3) 

V 2' 2 ; r(i)y4ifci?cosh 2 f V ; 

At large distances, we have cosh /j, ~ sinh /j^e M = ^ , the elliptic coordinate v is reduced to the ordinary polar angle 
(j>, and therefore the full asymptotic form of the solutions become 

* \/5 e ""~ ;> (c4 » 

From inspection of the large fj, behavior of the second solutions we find the asymptotic form to be 



Fey m (r) a y '-i-e^^+f) (C6) 



Gey m (r) « \j ~J {kr+ ^ (C7) 
Similarly we can also determine the asymptotic forms for the lower component radial solutions 



./,>,„(,) - \l ^~ im ^~ l[kr+l) (C8) 



Jr ,„(>-) v /^ e -^ e -i(*r+J) (C9) 



r f ,/,„„•) v /i e -*"- e -i(*r- f ) (C10) 
a, //,„(,) v /_l_ e -i»»- e -i(*r- f ) (C11) 



Appendix D: Plane wave expansion coefficients 



In this section we calculate the plane wave expansion coefficients and appearing in Eq. ( 94 ) 



* u ( ii ) - ist > {e-w-w ( ) ( ) } <»» 

In elliptical coordinates a unidirectional plane wave is described in terms of the incident direction 9 in the following 
way 

k = (cos Ox + sin 6y) , f = i?(cosh /i cos v x + sinh fj, sin v y) (D2) 
k-f — fci?(cosh jU cos f cos + sinh /i sin v sin 0) (D3) 
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Multiplying both sides by the transformation matrix (ST) 1 , 

^ ■ w ( z | ) - te<*-o ( ) + e ( ) 



fe. r - / cos § cosh f cos f + sin § sinh f sin § \ _ f Je^Be m \ f J 0r? Bo m 

6 { cos § cosh g sin | - sin § sinh f cos f ) ~ 1 2^ n m^ e ™W ^ . Je ^ Aem j + ^ i Jo^Ao m 



For brevity of notation, we write 



cos | cosh ^ cos | + sin | sinh ^ sin | \ _ / fupip-, v, 9) 
cos | cosh ^ sin | — sin | sinh ^ cos | 



We also have the following bi-orthogonal properties for the angular Whittaker Hill functions. 

27T 27T 



J Ao* s Be m du = k m S ms , J Ae* s Bo m du = k* m 8 ms (D5) 



o o 

27T 27T 



Ae*Be m dv = 0, J Ao*Bo m dv = (D6) 
o o 

where k m is the normalization constant. Operating on both sides with row vector ( Ao*(9) Ao* s (v) , Ao* s (9)Bo* s (v) ), 
integrating over v and 9 from to 27r, and applying the bi-orthogonality relations, we get 

J J gifcflCcoBhMooe^coBfl+BinhMBin^inflJ^.^ ( Ao;(l/)/ up (/i, I/, B 0* s (v) f d (^, V , 6) ) d V d9 = {< fc 2 ( Je s + * Je' s ) } 



Putting /n = m the above expression yields 

e _ / o 2T / o 2w e jfefl ( cosh ^ cos ^ cos0 + sinh ^ sini/sin8 ) J 4o*(6») ( Ao* s (v)f up ( l j l Q,u,6)+ Bo* s (v)f d (^,v,6) )dv d9 

^(JeM+zJe'M) (D7) 

Similarly, we can calculate n° by operating on both sides with row vector ( Ae* s (9) Ae* s (v) , Ae* s (9)Be* s (v) ) 

= J** C e l kR ( cosh "° cos v cos e+sinh ^' s ^^ 9 )Ael(9)(Ae* s ( V )f up (^, V ,9)+ Be* s ( V )f d (^, V ,9)) dvd9 
Hs ~ " " ~ " " _ " " (k* s r(JoM + Uo'M) [ ' 

We have checked that n% and n° are independent of the value of /xq. 
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